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Abstract 

The impact of a kination-dominated phase generated by a quintessential exponen- 
tial model on the thermal abundance of gravitinos and axinos is investigated. We 
find that their abundances become proportional to the transition temperature from the 
kination to the radiation era; since this temperature is significantly lower than the 
initial ("reheating") temperature, the abundances decrease with respect to their values 
in the standard cosmology. For values of the quintessential energy-density parameter 
close to its upper bound, on the eve of nucleosynthesis, we find the following: (i) 
for unstable gravitinos, the gravitino constraint is totally evaded; (ii) If the gravitino 
is stable, its thermal abundance is not sufficient to account for the cold dark matter 
of the universe; (iii) the thermal abundance of axinos can satisfy the cold dark matter 
constraint for values of the initial temperature well above those required in the standard 
cosmology. A novel calculation of the axino production rate by scatterings at low 
temperature is also presented. 
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1. Introduction 

A plethora of recent data [1, 2] indicates [3] that the two major components of the present universe are 
Cold Dark Matter (CDM) and Dark Energy (DE) with density parameters [1] 

(a) QcDM = 0.214 ± 0.027 and (b) Qde = 0.742 ± 0.03 (1.1) 



at 95% confidence level (c.l.). Identifying the nature of these two unknown substances, is one of the 
major challenges in contemporary cosmo-particle theories. 

The DE component can be explained by modifying the standard cosmology (SC) via the 
introduction of a slowly evolving scalar field called quintessence [4] (for reviews, see Ref. [5]). An 
open possibility in this scenario is the existence of an early kination dominated (KD) era [6], where the 
universe is dominated by the kinetic energy of the quintessence field; this period is an indispensable 
ingredient of quintessential inflationary scenaria [7, 8, 9]. During this era, the expansion rate of the 
universe is larger compared to its value during the usual radiation domination (RD) epoch. This implies 
that the reUc abundance of the weakly interacting massive particles (WlMPs) can be significantly 
enhanced with respect to (w.r.t) its value in the SC [10, 11, 12, 13], provided that they decouple from 
the thermal bath during the KD era. Further phenomenological implications of this effect for the future 
collider or astrophysics experiments have also been studied [14]. 

WIMPs are the most natural candidates [15] to account for the second major component of the 
present universe, the CDM. Among them, the most popular is the Ughtest neutralino [16, 17] which 



1 Introduction 



2 



turns out to be the lightest supersymmetric particle (LSP) in a sizeable fraction of the parameter space 
of sypersymmetric (SUSY) models and therefore, stable under the assumption of the conservation of 
/^-parity. However, SUSY theories predict the existence of even more weakly interacting massive 
particles, known as e-WIMPs [18] which can naturally play the role of LSP. These are the gravitino, 
G, and the axino, a (G is the spin- 3/2 fermionic SUSY partner of the graviton, and a the spin- 1/2 
fermionic SUSY partner of the axion which arises in SUSY extensions [19] of the Peccei-Quinn (PQ) 
solution [20] to the strong CP problem). As their name indicates, the interaction rates of gravitinos 
and axinos are extremely weak, since they are respectively suppressed by the reduced Planck scale, 
mp = Mp/ (Mp = 1.22 X 10^^ GeV being the Planck mass) and by the axion decay constant, 
fa ~ (IQi" - 10^2) GeV (for a review, see Ref. [21]). 

Due to the weakness of their interactions, e-WlMPs depart from chemical equilibrium very early 
(G at a energy scale close to mp and a close to fa) and we expect that their relic density (created 
due to this early decouphng) is diluted by the primordial inflation. However, they can be reproduced 
in the following ways: (i) in the thermal bath, through scatterings [22, 23, 24, 25, 26, 27] and decays 
[26, 28, 29] involving superpartners, and (ii) non-thermally [30, 31], from the out-of-equilibrium decay 
of the next-to-LSP (NLSP). In this paper we do not consider the possible non-thermal production of 
e-WIMPs, since this mechanism is highly model dependent (i.e., it is sensitive to the type and decay 
products of the NLSP). As a consequence, we do not consider either the out-of-equiUbrium decay of 
the one e-WIMP to the other (as in the case of Ref. [32], where G is the NLSP and a the LSP). In all 
these cases, extra restrictions have to be imposed in order not to jeopardize the success of the standard 
Big Bang Nucleosynthensis (NS). 

The latter requirement has to be satisfied also for unstable G. This restriction imposes a tight 
upper bound on the initial ("reheating") temperature, T\, of the universe in the SC [22, 23, 33, 34, 35]. 
On the other hand, if one of the e-WIMPs is a stable LSP, it has to obey the CDM constraint. In 
particular, its relic density Q-xh^ has to be confined in the region [1] 

(a) 0.097 < Q.xh^ < 0.12 for (b) 10 keV <mx< wnlsp (L2) 

where wtnlsp is the mass of the NLSP. Let us note, in passing, that the lower bound of Eq. (1.2a) is 
valid under the assumption that CDM is entirely composed by X's and the abundance of non-thermaly 
produced X's is negligible. The lower bound on mx arises from the fact that smaller mx cannot explain 
[36] the observed early reionization [1]. For about 10 < mx/^^eV < 100, X's may constitute warm dark 
matter (the mass limits above are to be considered only as indicative). 

In this paper we reconsider the creation of a KD era in the context of the exponential quintessential 
model [37, 38], taking into account restrictions arising from NS, the inflationary scale, the acceleration 
of the universe and the DE density parameter. Although this model does not possess a tracker-type 
solution [1 1, 39] in the allowed range of its parameters, it can produce a viable present-day cosmology 
in conjunction with the domination of an early KD era, for a reasonable region of initial conditions 
[40, 41, 42, 43]. We then investigate the impact of KD on the thermal production of e-WIMPs, 
solving the relevant equations both numerically and semianalyticaUy. We find that the abundance of 
e-WIMPs becomes proportional to the transition temperature from KD to RD era, Trr, and decreases 
w.r.t its value in the SC since Tkr can be much lower than Tj. In particular, we consider two cases, 
depending on whether T^r is higher - high T regime (HTR) - or lower - low T regime (LTR) - than 
a threshold Tc - 10 TeV, below which the thermal production of a via the decay of the superpartners 
becomes important; for this low temperature region, a novel formulae for the production of a's through 
scatterings is presented. It turns out that, in this part of the parameter space, a becomes an attractive 
CDM candidate within the quintessential kination scenario (QKS). 

Modifications to the thermal production of e-WIMPs have also been investigated in the context 
of extra dimensional theories [44, 45], where the expansion rate of the universe can be also enhanced 
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w.r.t its value in the SC, due to the presence of an extra term including the brane-tension. However, 
this increase is more drastic than in the QKS. In addition, the production of a via scatterings at low 
temperature and via the decay [26] of the SUSY particles has not been taken into account [45]. 

We start our analysis by reviewing the basic features of the exponential quintessential model in 
Sec. 2. We then present our numerical and semi-analytical calculations of the thermal abundance of 
e-WlMPs in Sec. 3 and study the parameter space allowed by several requirements for G (Sec. 4) and 
for a (Sec. 5). Our conclusions are summarized in Sec. 6. Computational issues on the low temperature 
a-production are discussed in Appendix A. 

Throughout the text, brackets are used by applying disjunctive correspondence, natural units 
(h = c = kB = 1) are assumed, the subscript or superscript refers to present-day values (except 
in the coefficient Vq) and log [In] stands for logarithm with basis 10 [e]. Moreover, we assume that the 
domain wall number [21] is equal to 1. 

2. The Quintessential Exponential IModel 

In this section we review the system of equations which governs the quintessential cosmological 
evolution (Sec. 2.1) and the various observational restrictions that we impose (Sec. 2.2). We then 
describe the salient features of this evolution in Sec. 2.3 and the allowed parameter space in Sec. 2.4. 

2.1. The Quintessential Set-up 

We assume the existence of a spatially homogeneous scalar field q (not to be confused with the 
deceleration parameter q in Sec. 2.2), which obeys the Klein-Gordon equation. In particular, 

q + 3Hq + V,q = 0, where V^Voe'^''"'' (2.3) 

is the adopted potential for the q field, the subscript , q [dot] stands for derivative w.r.t q [the cosmic 
time, t] and H is the Hubble expansion parameter, 

H = ^|Pq+pR+pM/ with pq = ]^q^ + V, (2.4) 

the energy density of q. The energy density of radiation, pR, can be evaluated as a function of the 
temperature, T, while the energy density of matter, p^, with reference to its present-day value: 

PR = ^gp* and pmR^ = P^i^o (2-5) 

(R being the scale factor of the universe). Assuming no entropy production due to domination of q or 
any other field, the entropy density, s, satisfies the following equations: 

sR^ = SpRl where s = -^9s* ■ (2.6) 

Here, the subscript "p" represents a specific reference point at which the quantities s and R are evaluated 
and gp*(T) [gs*(T)] is the energy [entropy] effective number of degrees of freedom at temperature 
T. The numerical values for these quantities are evaluated using the tables included in micrOMEGAs 
[49]. As it turns out, - g^t for T > Tns = 1 MeV, with values gp^ = 105.74 for T - 1 TeV 
and gpt, = 228.75 for T > 1 TeV assuming the particle content of the Minimal SUSY Standard 
Model (MSSM). Since the abundances under consideration take their present value on the eve of NS 
(T ~ Tns), we do not insist in a distinction between gfp* and Qs* and set gp* = gg* - g*. On the 
contrary, for T < Tns we obtain in general gp^ < gs*. Taking into account the existence of the cosmic 
background radiation at present, plus three(almost) massless neutrino species, we get [46] gp^, = 3.36 
and^°, = 3.91. 
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The numerical integration of Eq. (2.3) is facilitated by converting the time derivatives to 
derivatives w.r.t the logarithmic time [40, 43], which is defined as a function of the redshift z: 

T = In (/?//?())- -ln(l +z). (2.7) 
Changing the differentiation and introducing the following dimensionless quantities: 

Pm[R] = Pm[R]/Pc, = Vo/Pc and q = q/ V3mp, (2.8) 
Eq. (2.3) turns out to be equivalent to the system of two first-order equations: 

Q - Hq and HQ' + 3HQ + Vq = with ^p^+p^+p^, (2.9) 

where prime denotes (unless otherwise stated) derivative w.r.t x and the following quantities have been 
defined: 

V = Voe-^-*^, H = HI Ho, Q = Q/^, and p^ = + V. (2.10) 

In our numerical calculation, we use the following values: 

p° = 8.099 X 10"^^ /z^ GeV"* and Hq = 2.13 x W'^^h GeV (2.11) 

with h = 0.72. In addition, p^^ = 0.26 and Tq = 2.35 X 10'^^ GeV and from Eq. (2.5), we get 
p° = 8.04 X 10-^ 

Eq. (2.9) can be resolved numerically by specifying two initial conditions at a logarithmic time, Xj, 
which corresponds to a temperature Ti defined as the maximal T after the end of primordial inflation, 
assuming instantaneous reheating. We take ^(Xj) = throughout our analysis, without any lose of 
generality [13] and let as free parameter the value the square root of the kinetic-energy density of q at 
\ [13] 

y/p^ = Q(x,)/ - ^ - (2.12) 

The last equality holds with great accuracy, since we require a complete domination of kination at early 
times, as we describe below. 



2.2. Imposed Requirements 

We impose on our quintessential model a number of requirements which can be described as follows: 



2.2.1 . The Constraint of the Initial Domination of Kination. As we stress in the introduction, we focus our 
attention on the range of parameters that ensure an absolute or at least a relative initial domination of 
the ^-kinetic energy. This requirement can be quantified as follows: 

a\ = n^iTi) > 0.5 with = pg/ipq + PR + pm) (2. 13) 

the quintessential energy-density parameter. 



2.2.2. Nucleosynthesis Constraint The presence of pq has to preserve the successful predictions of Big 
Bang NS which starts at about x^^j = -22.5 corresponding to Tns = 1 MeV [47]. Taking into account 
the most up-to-date analysis of Ref. [47], we adopt a rather conservative upper bound on Q^(Xj,g), less 
restrictive than the one of Ref. [48]. In particular, we require: 

Qf - Q^(x^s) < 0.21 (95%c.l.) (2.14) 

where 0.21 corresponds to additional effective neutrinos species 6Ny < 1.6 [47]. In the left hand side 
of Eq. (2.14), we do not consider extra (potentially large [52]) contributions from the energy density 
of gravitational waves generated during a possible former transition from inflation to KD epoch [7]. 
The reason for this approach is that, inflation could be driven by another field different from q and 
therefore, any additional constraint from that period would be highly model dependent. Nevertheless, 
inflation may provide a useful constraint for the parameters of our model, as we discuss below. 
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2.2.3. Inflationary Constraint. Recent data [1] strongly favors the existence of an inflationary phase in 
the early universe. Assuming that this phase also generates the power spectrum of the curvature scalar 
Ps and tensor Pi perturbations, an upper bound on the inflationary potential Vi and consequently on 
Hi can be obtained, following the strategy of Ref. [53]. More specifically, imposing the conservative 
restriction r - Pt/Ps ^ 1, and using the observational normalization of P^ [1], we get 

Hi < —mpPll^ ^ Hi< 2.65 X 10^"^ GeV ^, < 1.72 x 10^^ (2.15) 
V2 

where * denotes that P^* is measured at the pivot scale = 0.002/Mpc. As we can see in Sec. 2.4 
the constraint on Hi is sufficient to restrict the parameter space of our model. Therefore, we are not 
obliged to impose other constraint on Tj, or equivalently on Ti, as was the case in Ref. [13]. 

2.2.4. Coincidence Constraint. The present value of p^, p^, must be compatible with the preferred 
range of Eq. (1.1b). This can be achieved by adjusting the value of Vq- Since, this value does not 
affect crucially our results (especially on the e-WIMPs abundances), we fix = p^/p^ to its central 
experimental value, demanding: 

Q°=p° = 0.74. (2.16) 

2.2.5. Acceleration Constraint A successful quintessential scenario has to account for the present-day 
acceleration of the universe, i.e. [1], 

-1 < WqiO) < -0.86 (95% c.l.) with Wq = (q'^/l - V)/iq'^/2 + V) (2.17) 

the barotropic index of the ^-field. In our case, we are not able to avoid [43] the eternal acceleration 
(u/g > -1/3 see below), which is disfavored by string theory. 

Let us also comment on our expectations for the redsift zt of the transition from deceleration to 
acceleration, and for the age of the universe, to. Due to observational uncertainties in the measurement 
of these quantities, we do not impose the data on them as absolute constraint. We estimate Zt by solving 
numerically the equation 

q(tt) = 0, with Zt = e""' - 1 and q = -RR/R^ ^-\-H'/H (2.18) 

the deceleration parameter; to is also estimated numerically, as follows: 

to= dR/R = (l/Ho) dx/Hix) (2.19) 

Jo Jxi 

We obtain 0.78 < Zt ^ 0.84 and 13.6 > ?o/Gyr > 13.4 as w^iO) varies in the range of Eq. (2.17). In 
both cases, as Wq(0) approaches - 1 the results on zt and to approach their values in the context of the 
standard power-law cosmological model with CDM and a cosmological constant (ACDM). The results 
on to are in agreement with observational data [1] (according to which ?o = (13.69 ± 0.26) Gyr at 95% 
C.L). On the other hand, our findings on zt are marginally consistent with the SNe la observations [2] 
according to which zt = 0.46 ± 0.26 at 95% c.l. However, it is probably premature to say more than 
that Zt is around unity (see, e.g., Ref. [54]). 

2.3. The Quintessential Evolution 

The quintessence field q undergoes three phases during its cosmological evolution [8, 13]: 
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Figure 1: The evolution of the quantities logp, with i = q (blacli lines), R+M (Ught gray line) and A 
(thin line) as a function ofT for A = 0.5, Ti = 10" GeV and Clf = 0.01 (^o = 2.04 x 10"). 



The kinetic-energy dominated phase during which pq is essentially given by q/2 » V; this results 
to u;^ ^ 1 and therefore, 

(where the last equality was extracted using Eq. (2.6)). The subscript or superscript "p" means that 
the various quantities are evaluated at a reference point p which can be the initial epoch (p = 1) or 
the epoch of the transition between KD and RD era (p = KR) or the eve of NS (p = NS). Applying 
Eq. (2.20) for p = NS, we derive the "transition" temperature Tkr from KD to RD era as follows: 

^NS\l/2^i _qNS\1/2 



P,iTKR)^PRiTKR) ^ TKR^T^s\^^j |^ ^^s j " ^^.21) 
For T > Tkr, the universe undergoes a KD epoch where H is adequately approximated by 

Applying this formula for p = NS for any given Ti, we get Hi as a function of Q^^. Therefore, 
we can use as a free parameter, instead of Hi. 

• The frozen-field dominated phase, where the universe becomes RD and pg is dominated initially 
by q/2 and subsequently by V. 

• The attractor dominated phase, where Pq - V dominates the evolution of the universe, and reaches 
the late-time attractor energy density: 

Pa ^ Pg (^/^o)"^^^*"'*^ with = - 1 (for ^ < VS) (2.23) 

the fixed point value of Wq. Today we obtain a transition from the frozen-field dominated phase 
to the attractor dominated phase [13, 40]. Although this is not a satisfactory solution to the 
coincidence problem, the observational data can be reproduced for a reasonable set of initial 
conditions, as we show below. 

The cosmological evolution described above is illustrated in Fig. 1 where we plot logp,- versus T 
for Ti = 10^ GeV (x, = -51.2), Q^s = o.Ol (H, = 7.07 x 10^^) and ^ = 0.5 (Vq = 2.04 x 10^^). For 
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log (T, / GeV) 

Figure 2: Allowed (gray and lightly grey shaded) region in the log Ti - log plane by Eqs. (2.13), 
(2.14) and (2.15). The conventions adopted for the various lines are also shown. 

i = q (bold black line), we show logp^, computed by inserting in Eq. (2.10) the numerical solution of 
Eq. (2.9). For / - A (thin black lines), we show logp^ derived from Eq. (2.23). For / = R + M (light 
gray line), we show log Pr+M' which is the logarithm of the sum of the contributions given by Eq. (2.5). 
For r > Tkr = 0.0098 GeV we have a KD era. We also obtain Wq(0) = 0.957 whereas = 0.916, 
Zt = 0.8 and to = 13.5 Gyr. 



2.4. The Allowed Parameter Space 

The free parameter space of our quintessential model is defined by the following quantities: 

A, Tj and or equivalently. A, T\ and Q^^. 

Agreement with Eq. (2.17) entails < i < 0.9 (compare with Ref. [43, 13], where less restrictive 
upper bounds on ^^(O) have been imposed)). This range is independent on t, and //, (or Ti and i^^^) 
as is shown in Ref. [13]. Given this fact, we focus on the last two free parameters of our model. 

The allowed area in the log T\ - log plane (which obviously is /l-independent), is indicated 
in Fig. 2. In the shaded areas Eqs. (2.13) (2.14) and (2.15) are fulfilled; in particular, in the dark [light] 
shaded area, = 1 [0.5 < < 1]. The left [right] boundary of the allowed regions (indicated by a 
dashed [soUd] line) is derived from Eq. (2.13) [Eq. (2.15)] whereas the upper boundary (indicated by a 
double dot-dashed line) comes from Eq. (2.14). We observe that the dashed and solid lines are almost 
straight. This can be understood by deriving the analytic relation between Ti and for fixed Hi or 
Q}q<l, namely: 



Ti - ^^fs 




and T, = Tms (2.24) 



"9 

The first [second] equation can be extracted by calculating Hi through Eq. (2.22) with p - NS [p = I 
and p = NS and eliminating Hi] and solving the resulting w.r.t T\. Taking the logarithms of the two 
sides of Eq. (2.24) and given that is quite small in the largest part of the available parameter space, 
we can convince ourselves that the dashed and solid curves in the log Tj-log plot have to be almost 
straight. Consequently, for a reasonable set of parameters {A, T\, fl^^), the exponential quintessential 
model can become consistent with the observational data, in agreement with [13, 40, 43]. 
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We assume that any relic abundance of an e-WIMP, X, due to its decoupling from the thermal bath, is 
diluted after inflation at an energy scale below 10^^ GeV. We compute its abundance produced through 

thermal scatterings and decays during the KD epoch and compare this result with the one obtained 
in the SC. The relevant equations are presented in Sec. 3.1. In Sec. 3.2 we describe our numerical 
evaluation and in Sec. 3.3 we derive useful approximate expressions. Their results are compared with 
the numerical ones in Sec. 3.4. 



3.1. The Boltzmann Equation 

The number density nx of X particles satisfies the Boltzmann equation, which can be collectively 
written in the following form [23, 26]: 

hx + 3Hnx = ^ J]{(v(r)ij}nf'nf + J]{ri}n^. (3.1) 

ij i 

Here, H is given by Eq. (2.4) and r^'^ is the equilibrium number density of the particle /, which can be 
adequately calculated in both the relativistic (m, <sc T) and non-relativistic (m, » T) regime [55] : 

«r - ^m]TK2{milT) (3.2) 

with Kn, the modified Bessel function of the second kind of order n and m,- [^,] the mass [the number 
of degrees of freedom] of the particle i. 

The quantity {{vcr)ij) is the thermal-averaged production rate of X from scatterings in the thermal 
plasma, the indices i and j run over all particles involved in the initial states of these processes, namely 

gluons ig), gluinos {g), quarks {q) and squarks {q), and the factor 1 /2 is introduced [23] in order to 
avoid double counting of the production processes. Formalistically speaking, the term {{vo-)ij)n^n^ 12 
corresponds to the so-called in Refs. [24, 25, 27] collision term Cx^f^^ and therefore, we can make the 
replacement: 

1 V N \ eq eq „ eq2 ^ / in the HTR, 

- ^<(.cr),,)„,. n. = Cxn with Cx = | ^ lTR, ^^'^^ 

and n*^ = ((3)T^/7t^, the equilibrium number density of the bosonic relativistic species (^(3) ^ 1.2 is 
the Riemann zeta function of 3). 

We proceed with the formulas for Cx's in the low (LTR) and high (HTR) T regime. In the 
relativistic regime (T » mi), Cx has been recently recalculated [24, 25, 27] in a consistent gauge- 
invariant treatment, using Hard Thermal Loop Approximation (HTLA) technics. The result is 



- 



(3.4) 



" ' in with^r^-— ^ for Z = a. 

Here, ga and Ma are the gauge coupling constants and gaugino masses respectively, associated with 
the gauge groups U{\)x, SU{2\ and 5f/(3)c, (fcj = (1.634, 1.312, 1.271) and (c„) = (33/5,27,72). 
Note that, contrary to the notation of Ref. [25], we use the Grand Unified Theory (GUT) inspired 
normalization of the hypercharge coupling constant g^. 

In our study, we calculate ga as a function of the temperature, solving the relevant one loop 
renormalization group equations with the MSSM particle content. In particular, we have 

-1/2 



1 



gaiT) = gcvT — ^ In — 



(3.5) 
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where (ba) = (33/5, 1, -3) and gcvT - 1/24 is the value of ^'s at the GUT scale Mqut - 2x 10^^ GeV. 
By adjusting ^gut and Mgut, we obtain the experimentally acceptable values of ga's at Mz- Below 
Mz, we use the gj running indicated in Ref. [56]. 

Throughout our analysis we impose universal initial conditions for the gaugino masses i.e., we 
assume that 

Ma(MGVT) = Ml/2, with a = 1,2,3. (3.6) 

The running of the gaugino masses Ma can be easily evaluated at one loop by solving the relevant 
renormaUzation group equations, which admit an exact solution: 



Ma(T) ^ MiJ-^] = Ml/2 (l - 1% In ■ 



(3.7) 



Mem I 

Under the assumption of Eq. (3.6), the lightest neutralino turns out to be a Bino (B), with mass 
nig ^ 0.41 Ml/2. 

From Eq. (3.4) one easily deduces that C^^ becomes negative for large enough values of g^. 

Indeed, from Eq. (3.5) we conclude that g^ increases as T decreases. In particular, C?^ > requires 

03 < 1.108 or T > Tc = 10^ GeV. Eq. (3.4), therefore, can be applied self-consistently for Ti > Tq 

and Tkr > Tq- Towards lower values of T, the finite masses m,, which have been neglected in deriving 

Eq. (3.4), and X-production from decays, start playing an important role. The correct inclusion of these 

effects is crucial in the case of a, since Q-ah^ takes cosmologically interesting values for T\ < Tq and/or 

Tkr < Tq, too. In Appendix A, we present a novel calculation of the term {(vcr)ij}nf^n^ in the case of 

J I J 

a and in the non-relativistic regime (T «: m,). Our final result can be cast in a form similar to this of 
Eq. (3.3), with: 



' y 

1675^(3)2 ^ 



ds Wij(s)Ki j pi(mi,mj) , 



(3.8) 



where - {mi + mj)^ and the symbols Wij and pi are defined in Eqs. (A2) and (A3). In particular, Wjj 
is related to the squared amplitudes of the various processes that contribute to {{vo-)ij), and are listed 
in Table 1. 

Finally, (F,) is the thermal-averaged rate of X-production from decays in the thermal plasma, 
which are related to the corresponding decay widths F,- through the formula [55]: 



<r,) = 



(3.9) 



KiimdT) 

These contributions have been recently calculated in the case of G in Ref. [29]; however, we do not 
include them in our calculation since their impact is roughly a factor of two. Such a minor change does 
not alter our results in any essential way. On the other hand, we do include these contributions in the 
case of a with i - g, q and B, applying the following formulas [28]: 



.9I9I 3 1 

3^«^3 



2\3 



for 



-ma 



fa 



mn In — 



m 



2 \ 



1-^ 



for i - q. 



for i = B. 



(3.10) 



50X32V/2 

In the above equations, g^, is calculated at a scale equal to m, if Ti > m,- or Ti if Ti < mi. Obviously, 
Tg and F^ are much more efficient than F^, due to the presence of g^. As we verify numerically, the 
contribution of these terms to the resulting D.xh^ becomes important [28, 26, 18] for T\ < Tsusy or 
^KR < TsusY, where Tsusy - mq. 
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3.2. Numerical Solution 

In order to find a precise numerical solution to our problem, we have to solve Eq. (3.1) together with 
Eq. (2.3). To this end, and following the strategy of Sec. 2.1, we introduce the dimensionless quantities: 

-nx = nxl (p?)'^' and rf^ = rf^l (p^)'^" . (3.11) 

In terms of these quantities (and substituting Eqs. (3.3) and (3.9) into Eq. (3.1)), this takes the following 
master form, for numerical manipulations: 

iin'^ + mux - Cxn"^^ - J] i'w^f Kx {nnlT) = 0, (3.12) 

where H is given in Eq. (2.9) and the following quantities have been defined: 

Cx = Cx V3mp (p«)'^' , m,- = f = and f = ^. (3.13) 

Eq. (3.12) can be solved numerically with the initial condition nx(Tj) ^ 0, where Xj corresponds to the 
initial temperature Tj. The integration of Eq. (3.12) runs from Tj down to Tns (an integration to also 
gives the same result). We find convenient to single out two cases: 

• In the HTR (Tj » Tq and Trr » Tq), we integrate Eq. (3.12) with Cx = Cf^ from Tj to Tns 
(in the integration for T < Tq vfQ take C^^ frozen at its value at T ^ Tq). As it turns out, 
contributions from Eqs. (3.8) and (3.10) are negligible. 

• In the LTR (Ti «; Tq or Tkr Tq), if (i) Ti > Tsusy. we integrate Eq. (3.12) successively 
from Ti to Tsusy with Cg = C?^ and then from Tsusy to Tns with Cg = whereas if 
(ii) Ti < TsusY we integrate Eq. (3.12) from Tj to Tns with = Ct^. It turns out that the 
contributions from Eq. (3.10) are important only for T < Tsusy- Note that for Tsusy < T < Tq 
(where the SUSY particles / with masses m,- ~ Tsusy are relativistic) there is no accurate result 
for {(v(r)ij), since neither HTLA is vahd (g^iT) > 1.1) nor our computation in Appendix A is 
apphcable (we consider non-relativistic particles i and j). However, we believe that the above 
procedure gives a result which is sufficiently accurate for our proposes, since small variation of 
TsusY (by 10%) leaves the result practically unaltered. 

The X yield, Yx - nx/s, and the relic density of X, 0.xh^ - p'j^h^/p^ (with px - mxnx), can be 
easily found, via the relations [57] : 

(a) yI = 3.55 X 10"^^ nx{T^s) and (b) ^xh^ = 2.748 x 10^ Y^ mx/GeV, (3.14) 
where [p^f^ /s(Tfis) = 3.55 x lO-^'^ and jq/p? = 2.748 x lOVOeV. 

3.3. Semi-Analytical Approach 

The crucial quantity Y^ for the computation of Q-xh^ via Eq. (3.14) can be also derived semi- 
anal ytically by employing a number of simplifications. We first re-express Eq. (3.1) in terms of the 
variable Yx = nx/s, in order to absorb the dilution term. Indeed [57] 

hx + 3Hnx = Yxs = Y'^Ts. (3.15) 
where prime in this section means derivative w.r.t T . Employing Eq. (2.6), we obtain: 

s = -3Hs =» r = -3Hs/s', (3.16) 
where H can be fairly approximated by applying Eq. (2.22) with p=KR: 
P^^^ I T \^ 
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The presence of > clearly indicates the deviation from the SC, where - and H = H^c = 
yfp^l VSmp. Inserting Eqs. (3.16) and(3.17) into Eq. (3.15) and substituting in Eq. (3.1), this can be 
rewritten in terms of the new variables as: 

Yx = Y^ + Yv (3.18) 

where 7^- and Yy satisfy the following differential equations: 

= -Mai ^Jl + ^g and Y^^-J] yY,KMT)IT' ^/TT^ . (3.19) 
In the above, we have defined the following quantities: 

jP- „eq 

Va = yx*Y^^Cx and yr, = ——^Y'^^yx*giTim^ with 7«l = — , (3.20) 

2^(3)^ s 

= 3^ = Vl-'^i? = ^ (' * ^) ^ 

The expression of as a function of Trr in Eq. (3.17) enables us to divide the integration of 
Eqs. (3.19) to two separated domains, since 1/yi + can be approximated as follows (compare with 
Refs. [44, 45]): 



^Jgf^/g* (Tkr/t) forr»rKR, 



i^^^-d forr«r;;: ^'-^^^ 

Following this simplification, Eq. (3.19) can be resolved trivially. Let us present our results for the 
HTR and the LTR separately. 

3.3.1 . The HTR. In this case {Ti » Tsusy and Tkr » Tsusy), we find that Y^ - 7^ where 

7« = P dT J^^y^ + f " rfr^^ (3.23) 

and y^^ - yo-{C^^). Fixing C^^ and g^ at their initial values in each domain of integration, an analytic 
formula can be derived, within an 10% accuracy: 

yo^iy^^^] in these, 

I ya'' ^g^/giTKR in (Ti/Tkr) + y^T^R intheQKS. ^ 
From the above expressions, we can easily deduce that Y^ decreases in the QKS w.r.t its value in the 
SC. This is to be expected, since in the SC it is proportional to Ti whereas, in the QKS, it becomes 
proportional to Tkr (lower than Ti). In the latter case, the first term (which is usually ignored in similar 
estimates [44, 45]) of the right hand side (r.h.s.) of Eq. (3.24), gives an equally important contribution. 

3.3.2. The LTR. In this case {Ti «; Tq or Tkr «: Tq), we find cosmologically interesting solutions 
only in the case of a. We focus on the most intriguing possibility, in which T\ » T^susy but 
7"kr ^susY- In this case, takes naturally a value close to its upper bound in Eq. (2.14), as 
can be verified via Eq. (2.21). We find 

yI = y^^^ + y^^^ + y^ where (3.25) 
I ■ (3.26) 

dTy'^J, (3.27) 
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with ]^ - i/o-(C^^). In the SC with T\ < Tsusy^ ^\ can be derived by summing the second 
terms of the r.h.s. of Eqs. (3.27) and (3.28) and setting Tkr = T\. Unfortunately, finding a general 
analytical expression for Eqs. (3.26) - (3.28) is not straightforward, mainly because the derivation of 
Ct^ from Eq. (3.8) requires a double integration of several lengthy squared ampUtudes (see Appendix 
A). However, for the benchmark values of m, used in our analysis: 

m~q = 1 TeV, = 1.5 TeV and mg = 0.3 TeV (3.29) 

we can write simple empirical relations which reproduce rather accurately our numerical results. We 
distinguish the following cases: 

• In the SC, the main contribution to 7? arises from the second term of the r.h.s. of Eq. (3.27) 
[Eq. (3.28)] with Tkr ^ T\ for > 0.3 TeV [Ti < 0.3 TeV]. Using fitting technics, we get a 
relation with a 15% accuracy: 

a-ah^ ^Am-ail+C Ti) e-'^l'^'lfl with A - 1.44 x 10^"* GeV, (3.30) 

B - 1A5A11 GeV and C = 0.001/GeV. The discrepancy can be attributed to the logarithmic 
factor involved in - see Eq. (3.10). This factor disturbs the dependence of on fa as 
written in Eq. (3.30). The origin of the exponential factor is the non-relativistic expansion [55] of 

^i( V^/r) {Ki{milT)\ involved in the definition of {{vcr)ij) [(r,)] - see Eq. (3.8) [Eq. (3.9)]. 

• In the QKS, the main contribution to F? comes from the first term of the r.h.s. of Eq. (3.27) 
and therefore, we expect that O-ah^ is independent of Ti. Our results can be reproduced from the 
following relation which holds with an excellent accuracy: 

a-ah^ ^ ^ = J^ma With D = 9.26 X 10^^ GeV (3.31) 

where we used Eq. (2.21) in the last equaUty. We observe that flg^^ oc 1/ .^Q^ or Qah^ oc Tkr. 

We conclude, therefore, that in both the QKS and the SC our calculation of is crucial in order to 
achieve a reliable result for Qg^^. 

3.4. Numerical Versus Semi-Analytical Results 

The validity of our semi-analytical approach can be tested by comparing its results for D.xh^ with those 
obtained by the numerical solution of Eq. (3.12). In addition, useful conclusions can be inferred for the 
behavior of Q.xh^ as a function of mx and £1^^. Our results are presented in Fig. 3. The lines are drawn 
applying our numerical code, whereas crosses are obtained by employing the formulas of Sec. 3.3. 
In Figs. 3-(a) [3-(b)], we display n^h^ [^ah^] versus [ma] for M1/2 = 0.7 TeV [fa = lO" GeV] 
and various O^^'s indicated on the curves. We take Ti - 10^ GeV, a value frequently met in the 
well-motivated models of SUSY hybrid inflation [58]. 

As we can verify via Eq. (2.21), the results for the HTR are applicable for any in Fig. 3-(a) 
and for Q^^ = 10'^^ in Fig. 3-(b), whereas for the residual Q^s in Fig. 3-(b) the results for UTH hold; 
in particular, the crosses are obtained from Eq. (3.23) [Eq. (3.25)] in the HTR [LTR]. In both regimes, 
Qxh^ decreases as increases, as we expect by combining Eqs. (3.23) and (2.21) for the HTR, 
and from Eq. (3.31) for the LTR. On the other hand, Q^/i^ cx l/m^, while ^5/1^ oc nta- This can be 
understood from the fact that 0.xh^ o<= mxY^ and 72 cx Cx or oc Ti. However oc \/m% (for 

G G 

Ma > m^) whereas C?^, and F,- are essentially independent on mg (see also Ref. [26]). 

In most cases, we observe that the semi-analytical findings approach rather successfully the 
numerical ones. Let us clarify, however, that the results presented in the following sections, are derived 
exclusively by our numerical program. 
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Figure 3: ilxh^ as a function of mx (X = G [X = a]) for various il^^'s, indicated on the curves, 
Ti = 10' GeV and Mi/j = 0.7 TeV [fa = 10" GeVJ (a /bj). ForCl^^ > lO"'', we take in our computation 
the values of ini indicated in Eq. (3.29). The solid lines [crosses] are obtained by our numerical code 
[semi-analytical expressions]. The CDM bounds of Eq. (1.2) are also, depicted by the two thin lines. 



4. Kination and Gravitino Thermal Abundance 

We first examine the impact of a KD era on the thermal abundance of G and study its cosmological 
consequences. We discriminate two cases, depending on whether G is unstable (Sec. 4.1) or stable 
(Sec. 4.2). 



4.1. Unstable Gravitino 

If G is unstable, it can decay after the onset of NS, affecting the primordial abundances of the light 
elements in an unacceptable way. In order to avoid spoihng the success of Big Bang NS, an upper 
bound on 7^ is to be extracted as a function of and the hadronic branching ratio of G, Bh 
[33, 34, 35]. Let us specify some representative values of this constraint, taking into account the 
most up-to-date analysis of Ref. [33]. In particular, if G decays mainly to photon and photino, from 
Fig. 1 of Ref. [33] we can deduce: 

p60GeV 

Y^iTm) ^ \ 10"^^ for ^ \ 600 GeV and Sh = 0.001, (4.1) 
[ 10~i3 [ 11 TeV 

whereas if G decays mainly to gluons and gluinos, from Fig. 2 of Ref. [33] we can deduce: 

( 10-1^ poo GeV 

J'g(^ns) ^ ] 10"^^ for ^ \ 680 GeV and Bh = 1. (4.2) 
[9.6x10-15 [ 10 TeV 

We observe that for Sh - 1, the upper bound on Y^{Tt^s) does not exclusively increase with an increase 
of Mq, as in the case for Bh = 0.001. 

In the SC (where no late-time entropy production is expected), setting M1/2 = 500 GeV, we obtain 
a stringent upper bound on Tf. 

( 2.3 X 10^ GeV ( 360 GeV 

Ti < i 4 X lO'^ GeV for ^ \ 600 GeV and Bh = 0.001, or (4.3) 
[6x10^ GeV [ 11 TeV 

['8.5xl05GeV pOO GeV 

Ti < i 3.1 X 10^ GeV for ^ i 680 GeV and Bh = 1. (4.4) 
[5.4x10^ GeV [ 10 TeV 
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Figure 4: Areas in the logTi .^^—^ 
requirements (gray and lightly gray shaded area) 



log Q^^ plane that are simultaneously allowed by the quintessential 



Eqs. (2.13), (2.14) and (2.15) - and the gravitino 



constraint for = M1/2 = 500 GeV and Sh = 0.001 (black lined area) or 8^ = 1 (white lined area). 



Clearly the upper bound on Ti becomes significantly more restrictive for large B^s and low m^'s. 

In the QKS, if we set Ti = 10^ GeV and M1/2 = 500 GeV we can obtain a lower bound on Q^^ 
which can be transformed to an upper bound on Tkr via Eq. (2.21). In particular, 



ForSh = 1, 



0.001, 








6 X 10" 


-19 




1 X 10" 


-21 




.6x 10- 


-25 


1, 







5 X 10-1^ 
5.5 X 10-1^ 
, 4.5 X 10-22 



2.8 X 10^ GeV 

6.8 X 10^ GeV for - 

2.8 X 10^ GeV 



9.6 X 10^ GeV 
4.3 X 10"^ GeV for 
m'^ GeV 



360 GeV, 
600 GeV, 
11 TeV. 

200 GeV, 
680 GeV, 
10 TeV. 



(4.5) 



(4.6) 



The importance of a KD era in avoiding the gravitino constraint can also be induced by Fig. 4, 

where, in contrast to our previous approach, Ti is now variable, whereas is fixed to a representative 
value. In Fig. 4, we show the regions in the log - log plane that are allowed by the quintessential 
requirements (see Fig. 2), for = M1/2 = 500 GeV and Bh = 0.001 (black lined area) or 5h = 1 
(white lined area). We observe that for Sh = 0.001 the required minimal is lower that in the case 
of Sh = 1. This is because for Sh = 0.001 we impose Yq(Tj^s) ^ 2 x 10"'^, whereas for Sh = 1, 
we impose Y^{Tj<is) ^ 2 x 10""'^ (in accordance with Figs 1 and 2 of Ref. [33]). As a consequence, 
the maximal allowed T^r ^ (4.8 x 10^ - 6.8 x 10^) GeV for Sh = 0.001 is higher than the one 
(4.6 x 10"^ - 5.7 x 10^) GeV allowed for Sh = 1. Finally, we observe that the minimal depends 
very weakly on Ti. 

We can conclude, therefore, that the gravitino constraint can be totally eluded in the QKS, even 
for tiny values of D^^. 
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Figure 5; Allowed (lightly gray shaded) regions in thew-Q-Ti plane for G-CDMwithO.5 < Mi/2/TeV < 
1 in the (a) SC (Q^^ = 0) and (b) QKS with = lO'^^ The conventions adopted for the various lines 
are also shown. 



4.2. Stable Gravitino 

G can be stable if it is the LSP within SUSY models with R-parity conservation. In this case, 
it constitutes a good CDM candidate, provided D.^h^ satisfies Eq. (1.2). Imposing this condition 
constrains the free parameters, which in this case are mg, Ti, M1/2 and Q^^. 

In Fig. 5 we present the allowed regions (lightly gray shaded) in the - T\ plane in the (a) SC 
(Q^^ = 0) and (b) QKS with = lO'^^ letting M1/2 vary in the interval (0.5 - 1) TeV. The black 
solid [dashed] lines correspond to the upper [lower] bound on Q.^h^ in Eq. (1.2a), whereas the gray 
solid lines have been obtained by fixing n^/i^ to its central value in Eq. (1. 2a) for M1/2 = 0.7 TeV. In 
practice, the soUd [dashed] line is constructed for M1/2 = 0.5 TeV [M1/2 = 1 TeV]. This is because 
Q-^h^ oc Ti Mllm^ [D.^h^ oc In Ti Mllm^l in the SC [QKS] as deduced by Eqs. (3.24) and (3.4). 

The upper boundary curve (dot-dashed line) in Fig. 5-(a) results from the requirement that G is 
the LSP and thus is lighter than B. Note that varies from 205 to 411 GeV for M1/2 in the interval 
(0.5 - I) TeV. On the other hand, the upper boundary curve (double dot-dashed Une) in Fig. 5-(b) 
comes from the upper bound on in Eq. (2.15). The lower boundary curve (dotted [thin] line) in 
Fig. 5-(a) [Fig. 5-(b)] arises from the saturation of ^73 < 1 [Eq. (2.13)]. Note that towards lower T's 
the preferred tends to almost unnaturally small values. We conclude therefore, that the LTR is not 
cosmologically interesting in the case of G-CDM. 

From Fig. 5-(a) we observe that in the SC (O^^ = 0) somehow larger Tj's than in the case with 
unstable G are allowed (compare, e.g., with Eq. (4.3)). Much larger Tj's and m^'s are allowed in the 



QKS - see Fig. 5-(b). However, in this case, we are obliged to fine tune to a very low value 10 
(^KR = 7 X 10^ GeV), in order to obtain acceptable Q^/z^. We consider such a low as unnatural 
and thus, conclude that G is not a good CDM candidate within the QKS. 
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5. KiNATION AND AxiNO THERMAL ABUNDANCE 

We now turn to a, as the main candidate for CDM in the universe. The production of a is usually 
accompanied by the production of a scalar SUSY partner, the s-axion; this may undergo out-of- 
equilibrium decays [59], producing entropy and therefore, diluting any prexisting Qq/i^. In our 
analysis, we assume that this is not the case [26]. Moreover, we limit ourselves to Ti < fa since 
for larger Tj's the PQ symmetry [20] is restored and so, no particle from the axion supermultiplet has 
been produced. 

Let us initially derive the temperature To at which a decouples from the thermal bath in the SC 
and QKS. The a-decouphng occurs when 

^(^d) - raCTo), where = {(vcrW ~ 6 A^F(A^f - l)gigj (5.1) 

is an update [60] of the interaction rate which maintains a's in chemical equilibrium with the thermal 
bath (we take A'^p - 12 and A^3 = 3 as explained in Appendix A). Similarly to the previous discussion, 
we distinguish two cases: 

• In the SC, solving Eq. (5.1) w.r.t (after replacing H by its expression in Eq. (3.17) for = 0), 
we find 

6 X 10^ ^ r^^/GeV ^ 2 X lO^'^ for 10^° < fa/GeV < 10^^. (5.3) 

• In the QKS (after replacing H in Eq. (5.1) by its expression in Eq. (3.17) with given by 
Eq. (3.22)), we find 

H > r, for Tkr < or > ^ f^]' (5.4) 

C 1 X 10-^^ for fa = 10^° GeV, 

[ 1 X 10-28 for fa = IO12 GeV. 
These values were extracted numerically without the approximation of Eq. (3.22) and therefore, 
are less restrictive that the result of the analytical estimate in Eq. (5.4). Note that the constraints 
on are independent on Tj. 

As we already emphasize in Sec. 3 we expect that the primordial (i.e., due to the a-decoupUng) a 
yield [27] (yP = (3n^^/2s)(TD) - 1.8 x 10"^) is diluted by the entropy release during reheating of the 
universe to a temperature Ti. Under this assumption, we have to ensure that Ti < [Tkr < T^] 
constructing the regions where a is a viable CDM-candidate in the SC [QKS]. 

Our results for the HTR and the UTR are analyzed separately in the following. Let us remind that 
the discrimination between the two regimes (HTR or LTR) depends on the hierarchy not only between 
Ti and Tc but also between Tkr and Tc- In the SC we have Ti » Tq [Ti «; Tc] for the HTR [LTR] 
(see Sec. 5.1 [Sec. 5.1]). In the QKS we take Ti » Tq and Tkr > Tq for the HTR (see Sec. 5.1), but 
we consider as more natural choice (motivated by the majority [58] of the inflationary models) to take 
Ti » Tq and Tkr < Tq for LTR (see Sec. 5.2). The dependence of our results on the variation of Tj in 
the QKS is also studied in Sec. 5.2. 

5.1. The High T Regime 

In this case the HTLA works well and therefore the comparison with the case of G-CDM is 
straightforward; the free parameters in the present case are: mg, fa, Ti and D^^. 
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Figure 6: Allowed (lightly gray shaded) regions for a-CDM in the HTR and nta — Ti plane with 
10'° < /a/GeV < 10'^ for the (a) SC (flf = 0) and (b) QKS with Cl^^ = lO"". The conventions 
adopted for the various lines are also shown. 

In Fig. 6 we display the allowed regions (lightly gray shaded) in the ma - Ti plane in the (a) 
SC (O^s = 0) and (b) QKS with Q^s = iQ-is, letting vary in the interval (10i° - lO^^) GeV. The 
selected is the largest possible value that ensures the validity of HTLA, since Tkr ^ 10^ GeV. The 
black solid [dashed] lines correspond to the upper [lower] bound on Oa/i^ in Eq. (1.2a), whereas the 
gray solid lines have been obtained by fixing Qah^ to its central value in Eq. (1.2a) for fa - 10^^ GeV. 
In practice, the solid [dashed] line is constructed for /„ = 10^^ GeV [fa - lO'** GeV]. This is because 
n^h^ oc Ti m-alfl [Q.ah^ oc In Ti nia/f^] in the SC [QKS] as deduced by Eqs. (3.24) and (3.4). 

The left boundary curve (dot-dashed line) in Fig. 6-(a) comes from the lower bound of Eq. (1.2b). 
On the other hand, the upper boundary curve (double dot-dashed Une) in Fig. 6-(b) comes from the 
upper bound on Ti, Tj < fa. The lower boundary curves (dotted lines) in both Fig. 6-(a) and Fig. 6-(b) 
come from the saturation of 93 > 1 which determines the range of validity of the HTLA. Needless to 
say that, for Ti < 10"* GeV, imposing Eq. (2.13) also fails, since Ti - Tkr. 

We observe that in the SC the lower bound of Eq. (1.2b) is more restrictive than our requirement 
Ti < T^, which is satisfied for all fa's - see Eqs. (5.2) and (5.3). On the other hand, in the QKS, the 
selected Q^^ fulfills Eq. (5.5). For central values of and Qah^ we find 1.85 x 10^ < Ti/GeV < 10"^ 
for 10-5 < ma/GeV < 0.018 in the SC, whilst 10^^ < Ti/GeV < 10"^ for 0.0015 < m^/GeV < 0.02 in 
the QKS. We observe that in the SC the allowed Ti and are rather low, whereas in the QKS larger 
Ti and mg are permitted. Note, also, that towards lower Ti the preferred tends to larger values, 
implying that the LTR is cosmologically interesting in the case of a-CDM (contrary to the case of 
G-CDM). 

5.2. The Low T Regime 

In this regime Q-ah^ includes sizeable contributions from thermal scatterings and decay of SUSY 
partners. As a consequence, our computation depends not only on the parameters of the previous 
case (mg, fa, Tj and n^^), but also on the masses of the SUSY particles, m,. We use the m,'s of 
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Figure 7: Allowed (lightly gray shaded) regions for a-CDM, 10'" < /a/GeV < 10'^ and mi 's shown in 
Eq. (3.29) in the LTR and (a) in the m-a - Ti plane for the SC (fl^^ = 0) and (b) in the m-a - Q^^ plane for 
the QKS with Ti = 10' GeV. The conventions adopted for the various lines are also shown. 



Eq. (3.29) for both the SC and the QKS. Needless to say that our prehminary assumptions Ti < fa and 
Ti < r^c [Tkr < r^C] are comfortably satisfied in the SC [QKS]. 

In the SC (H^^ = 0), the allowed (lightly gray shaded) region for a-CDM is shown in the nia - Ti 
plane - see Fig. 7-(a). On the other hand, in the QKS, our results are independent of Ti - see Eq. (3.31); 
consequently, we depict the allowed (lightly gray shaded) region in the - Q^^ plane with fixed 
Ti = 10^ GeV - see Fig. 7-(b). The black solid [dashed] line corresponds to the upper [lower] bound on 
Qah^ in Eq. (1.2a), whereas the gray solid lines have been obtained by fixing Qah^ to its central value 
in Eq. (1.2a) for = 10^^ GeV. In practice, the solid [dashed] fine is constructed for fa = 10^^ GeV 

[/^ ^ 10"' GeV]. This is because Ti oc Blnif^Qah^) - Sln(Ama) [-^[^ oc ma/nah^f^] in the SC 
[QKS] as deduced by Eq. (3.30) [Eq. (3.31)]. 

In the SC, the lower bound of Eq. (1.2b) determines a part (dot-dashed line) of the left boundary 
curve in Fig. 6-(a), whereas the upper bound of Eq. (1.2b) sets an upper bound (thin fines) in both 
Fig. 6-(a) and Fig. 6-(b). The upper boundary curve (double dot-dashed line) in Fig. 6-(b) comes from 
the upper bound on in Eq. (2.14). The allowed area in Fig. 6-(a) [Fig. 6-(b)] terminates from above 
[below] at Ti - 1 TeV [T^r - 1 TeV], so that our formulas for Ct^ to be fully applicable. 

A sharp suppression of ^a^^ is observed, due to the Boltzmann suppression factor e~^l^^, in the 
SC - see Fig. 7-(a). In this case, larger but very low T\ are allowed. On the contrary, T\ can be 
fixed to a naturally [58] high value within the QKS. At the same time, and fl^^ take interestingly 
large values in the allowed region - see Fig. 7-(b). For these reasons, this case is considered as the 
most intriguing of this paper. For central values of fa and Qg/z^ we find 58.5 < Ti/GeV < 1000 for 
300 ^ ms/GeV ^ 8 x 10"^ in the SC, whereas lO-^^ ^ Q^s ^ o.056 for 8.5 x 10"^ ^ mg/GeV ^ 300 
in the QKS {Ti = 10^ GeV). Our formalism of Appendix A is of crucial importance in order to obtain 
a reliable result in both the SC and QKS. 

Finally, it would be interesting to directly compare the naturality of the a and G as CDM 
candidates in the QKS. This is done in Fig. 8 where, contrary to our strategy in Fig. 5-(b) [Fig. 7- 
(b)], we fix [mg] to some exemplary value and let Ti and vary in their allowed region of 
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log (T, / GeV) 



Figure 8: Simultaneously allowed areas in the log Ti - log fl^^ plane taking into account the 
quintessential requirements (gray and lightly gray shaded area) - Eqs. (2.13), (2.14) and (2.15) - and 
the CDM constraint for G-CDM (black lined area) with = 100 GeV and 0.5 < Mi/i/TeV < 1 or 
a-CDM (white lined area) with m-a = 5 GeV, m, 's ofEq. (3.29) and lO'" < /JGeV < 10'^ (white lined 
area). The bold black [white] line is obtained for M1/2 = 0.7 TeV [fa = 10" GeVJ and corresponds to 
Qxh^ - 0.11 (with X = G[X = a])- 



Fig. 2. In particular, in Fig. 8 we present in the log Ti - log D.^ plane the allowed region by both the 
quintessential requirements and the CDM constraint for G-CDM or a-CDM. The first set of constraints 
- Eqs. (2.13)-(2.15) - is satisfied in the gray and lightly gray shaded area. The CDM constraint - 
Eq. (1.2) - for G-CDM with = 100 GeV and 0.5 < Mi/2/TeV < 1 is fulfilled in the black 
fined region, where Q^^ is tuned to low values. The upper [lower] boundary curve of the black 
lined region corresponds to Q^h^ - 0.097 [^^h^ = 0.12] and is constructed for M1/2 = 1 TeV 

[Ml/2 = 0.5 TeV]. This is because .^j^^ oc Ml/m^Q^h^ as concluded by Eqs. (3.4) and (3.24). 

On the contrary, Eq. (1.2) for a-CDM with mg = 5 GeV and 10^° < fJGeV < lO^^ is met in the 
white lined region, with much more natural O^^'s than the ones needed in the black lined region. 

For the required Q^^'s our results of the LTR are applicable. For Ti > I TeV, Q^/i^ is obviously 
independent on Ti and Eq. (3.31) approaches fairly our numerical results. The upper [lower] boundary 
curve of the white lined region corresponds to D.ah^ = 0.097 [Og/i^ - 0.12] and is constructed for 
fa - 10^^ GeV [fa - 10^^ GeV] (in accordance with our discussion above). The bold black [white] 
fines are constructed for M1/2 = 0.7 TeV [fa = 10^^ GeV] and correspond to Q.xh^ - 0.11. For these 
values, we find 8 x 10*^ < Tj/GeV < 2 x 10" and 6.3 x lO'^'' < Q^s < ^9-23 for G-CDM, whereas 
92.2 < Ti/GeV < 1.8 x 10'^ and 1.25 x lO^'^ < Q^s < 1.6 x lO'^ for a-CDM. Therefore, in the 
framework of the QKS, a is clearly a more natural CDM candidate than G. 

6. Conclusions 

We presented an exponential quintessential model which generates a period dominated by the kinetic 
energy of the quintessence field. The parameters of the quintessential model (A, Ti, Q^^) were confined 
so that 0.5 < 0.q(T{) < 1, and were constrained by current observational data originating from NS, the 
acceleration of the universe, the inflationary scale and the DE density parameter. We found < A < 0.9 
and studied the allowed region in the (Ti, Q^^)-plane. 

We proceeded to examine the impact of this KD epoch to the thermal abundance of G and a. We 
solved the problem (i) numerically, integrating the relevant system of the difi'erential equations equation 
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that governs the evolution of the X-number density and (ii) semi-analytically, producing approximate 
relations for the cosmological evolution before and after the transition from KD to RD, and solving 
the appropriately re-formulated Boltzmann. Although we did not succeed to achieve general analytical 
solutions in all cases, we consider as a significant development the derivation of a result by solving 
numerically just one equation, instead of the whole system above. Moreover, for typical values of m,'s 
in Eq. (3.29), empirical formulas that reproduce quite successfully our numerical results were derived. 
For unstable G, the G-constraint poses a lower bound on Q^^, which turns out to be almost independent 
of Ti. The CDM constraint can be satisfied by the thermal abundance of G for extremely low values 
(10~^^ - 10"^"*) of Q^^. On the contrary, the former constraint can be fulfilled by the thermal abundance 
of a with values of Cl^^ close to the upper bound posed by the requirement for successful NS. 

Let us also comment here on three minor subtleties of our calculation which, do not alter the basic 
features of our conclusions (although could potentially create some quantitative modifications to our 
results). In particular: 

• Throughout our investigation we did not identify the nature of NLSP. Therefore the upper bound, 
shown in Fig. 5-(a) [Fig. 7], on [mg] derived from the requirement [m^ < ntg] nia < could 
be modified if there is another SUSY particle fighter than B. Moreover, we did not consider the 
NS constraints concerning the late decays of the NLSP into X\. These constraints [26, 30, 31] 
depend very much on the properties of the NLSP, i.e. its composition, its mass relative to the 
mx, and its coupling to Z's. Consequently, additional bounds on the mx might arise. As the a 
interactions are not as strongly suppressed as the G interactions, the a LSP anyhow is far less 
problematic than the G LSP w.r.t these constraints. 

• In the case of G, we did not incorporate contributions to Q^/i^ from the process of reheating. 
Indeed, these extra contributions can be a fraction of the result shown in Eq. (3.24) [34, 62], in 
the case of the usual reheating realized by the coherent oscillations of a massive particle [61]. 
However, in the QKS several reheating processes have been proposed [63] and therefore, any safe 
comparison between the SC and the QKS has to be performed for T < Tj. In other words, to 
keep our investigation as general as possible, we preferred to study the evolution of the universe 
after the start of the RD [KD] era in the SC [QKS] (we simply assumed the existence of an earlier 
inflationary epoch). 

• In the case of the a-CDM, we used throughout our investigation some representative masses for 
the superpartners, given in Eq. (3.29). Variation in these values (especially in m§ and m§) has 

an impact on Og/i^ in the LTR (e.g., for nig - tUq the contribution of the process q*q to can 
be enhanced). In addition, a further uncertainty in our calculation arises from the determination 
of TsusY below which C?^ is replaced by in the integration of the relevant equations. Note 
however, that the aim of the present paper is to demonstrate the change of O-xh^ due to the presence 
of the KD era and not a full scan of the SUSY parameter space. 

Although our results have been derived in the context of an exponential quintessential model, their 
appficabifity can be extended to every model that generates a KD phase, even without [6] quintessential 
consequences. It is worth mentioning that in the presence of kination we can obtain simultaneous 
compatibility of both the gravitino and the CDM constraint (i.e., the lower bound on Q^^ from the 
gravitino constraint is compatible with the Q^^'s needed in order to have the correct amount of a 
CDM). It would be probably interesting to check if we can obtain a simultaneous compatibility of these 
two constraints with additional bounds, arising e.g. from leptogenesis and neutrino masses [64, 65] or 
the quintessino abundance [66]. 
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A. AxiNO Production via Scatterings at Low Temperature 

In this appendix we present the necessary ingredients used for the evaluation of the a-production via 
scatterings at low T . Our starting point is the general formula in Ref. [67], which gives the thermal 
averaged cross section times the relativistic invariant relative velocity, {{vo-)ij}, of two particles / and j 
with masses m,- and mj (in general m,- ^ mj) and degrees of freedom gi and gj respectively: 



where = (w; + mjf and 



"''^ 32^ Ji ^'^"^^cmIMvI • (A2) 

Here, \Mij\^ is the squared amplitude of the process / + j ^ k + a which contributes to {{vo-)ij), nik 
is the mass of the outgoing particle k in this process, ^cm is the scattering angle in the center-of-mass 
frame and pi and pf are the magnitudes of the incoming and outcoming 3-momentums in the same 
reference frame; these are given respectively by: 

Kmumj) Umk,ma) 

Pi(mi,mj) — and pf(mk,ma) p — (A3) 

ys ys 

where X^(m,-, m^) -^{s- ("^i + m^)^) - (m, - m^t)^) (A4) 

Multiplying {(vcr)ij} by n^^nj^ and using Eq. (3.2) we arrive at our final result in Eq. (3.3), with 
Cx^CLT given by Eq. (3.8). 

For the manipulation of this result it would be useful to remember that \Mij\^ is in general a 
function of the Mandelstam variables, s, t and u (and the masses of the involved particles). However, 
one of the Mandelstam variables, typically u, can be ehminated in favor of the other two, through the 
formula: 

u - + mj + mj^ + - s - t (A5) 

whereas in the reference frame of the center of mass, t can be expressed in terms of s and 6cm using 
Eq. (A5) and the relation [67]: 

(mf - mj)(ml - mj) 

t-u- 1- Api{s)pi{s) cos ^CM (A6) 

s 

Therefore, lAl^p can be written as a function of only s and ^cm- 

The computation of the various |Al,yP can be realized applying standard techniques [68]. We 
concentrate on the a production processes which involve SU{3)c interactions. The Feynman rules for 
these interactions are indicated in Fig. 9 and originate from the following Lagrangian term [26, 27]: 

2 

£~a~gg = « y 5 75 {f ,f] f G^^, where ga = ' (A7) 
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'abc 



Figure 9: Feynman rules used in tlie derivation of\Mij\~ for the vertices wifh a. Tiie arrows were drawn 
according to the conventions of Ref. [68]. 

^fiv - (^M^v ~ ^v^p - Oif^^^g^d'^fi is the gluon field tensor, /^""^ are the structure constants of the SU{3)c 
algebra with a, b, c = 1, 8 whereas fi and v are the usual spacetime indices. The residual Feynman 
rules needed for our calculation are taken from Ref. [69]. In Table 1, we list the Feynman diagrams 
included in the calculation of |At,yp. Note that in the non-relativistic regime some of the processes {gg, 
qq and gq) contributing to lAl/yP in the relativistic regime [27] are kinematically blocked and thus, are 
non included in our computation. 

In Table 2 we list analytical expressions for the non-averaged squared amplitudes gigj\Mij\^'s of 
the various scatterings of ij particles. Let us make some comments on the findings that we present in 
this table: 

• Although our numerical results for the mfs of Eq. (3.29) are practically independent on ma (since 
it is expected to be significantly lower than ni,'s - especially rUg and we explicitly display it 
in our formulas. Only gluons and quarks are taken to be massless. Due to the smallness of ma 
as compared to /m,'s, Ct^ can be approximated by an interpolated function of T, for every chosen 
m/'s. This is, in practice, a great simplification for our numerical treatment. 

• The Feynman gauge has been used throughout our computation. The substraction of the 
unphysical mode in the process gg with two external g's can be obtained by employing appropriate 
projection operators [70] for the transverse polarization states. 

• We present the complete expressions for all gigj\Mij\^'s, except for the case of gg, where 
gigjlMijl^ turns out to be very lengthy; in this latter case, we display only the important 
contributions which arise from the t-channel g exchange and the interference terms between the t 
and s channel and PI. We have checked that our results on the total squared matrix elements for a 
production agree nicely with those of Ref. [27] in the relativistic limit (setting mg = mq - ma - 0). 

• In the expressions we present, ^^^ylAI/yP's are weighted with appropriate multiplicities following 
Ref. [69]. In particular, in the second column of Table 2 we arrange factors due to the summation 
over all the color-triplets {N-p = 12) - taking into account the corresponding charge conjugate 
multiplets which imply an extra factor of 2. The factor 1 /2 occurs in the gg process because of 
the identical particles in the initial state. Factors arising from the summation over color degrees 
of freedom and color indices, according to the identities 

|/''f = ^3(A^3' - 1) and 1^// = \ (^3 - l) . (A8) 

a,b,c a,I,J 

are included clearly in the expressions of the third column {T^j are the generators of the SU{3)c 
group and A^3 - 3). 

• The decay width of g is evaluated employing the tree level result of Ref. [71] - possible next-to- 
leading order corrections are subdominant. The result is checked to ensure that it is consistent 
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Final 


Interaction 


State 
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sig""), t{g\ u{g^) 




qj a 


Kg"") 


g'^-g" 


f a 


s{g% Kg''), u{g% PI 


qi9^ 


qja 


tig"") 


qiq*j 


~g^~a 


s(g') 


qi g^ 


qj a 




q*i qj 




sir) 



Table 1: The Feynman diagrams contributing to {ivcr)ij) with a in the Snal state. The exchanged particles 
are indicated for each relevant pair of initial and final states. The symbols s{i), t{i) and u{i) denote tree- 
graphs in which the particle i is exchanged in the s-, t- or u-channel and PI stands for "point interaction". 
The superscripts are S U{3) color indices, whereas the subscript I and J are family indices. 



with the output of the calcHEP package [72], and is the following: 



(A9) 



Although some integrals in Eq. (3.8) (especially the ones corresponding to the processes gg and 
qg) turn out to converge slowly, we checked that our results are pretty stable without the need to 
introduce an effective mass for the gluons propagators, as in the relativistic case [24, 26, 27]. 

In our numerical computation and with the m,'s listed in Eq. (3.29), we do not include the 
contributions of the processes qq* and qg since they are, in general, negligible. The major 
contributions to |AI,7p come from gg and qg for T - Tsusy and from qq for T «; Tsusy- The 
processes gg and qg give in general important contributions which cannot be neglected. Note, 
however, that the contribution of the process q*q could also be enhanced, if ntg ^ ntq, since we 
would have a resonance in the ^ - ^ annihilation via an j-channel exchange of a g. However, we 
consider that the choice nig - niq would be a rather ugly tuning of our free parameters and thus, 
we opted to use the m,'s in Eq. (3.29) in order to demonstrate our findings. 
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Table 2: The quantities gigjlMijl^ for a production from scatterings of the ij particles in the non- 
relativistic regime. The results are summed over spins, generation and color indices, in the initial and 
final state. The Mandelstam variables are defined as s = {p\ + p2)^, t = (pi - kiY and u = (pi - ^2)^ 
where the four-momenta pi, p2, ki, and k2 are associated with the particles in the order in which they 
are written down in the columns "initial state" and "final state" of Table 1. 
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